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Abstract. The Kibble-Zurek mechanism is applied to the spontaneous 
formation of vortices in a harmonically trapped thermal gas following a 
temperature quench through the critical value for Bose-Einstein condensation. 
While in the homogeneous scenario vortex nucleation is always expected, we show 
that it can be completely suppressed in the presence of the confinement potential, 
whenever the speed of the spatial front undergoing condensation is lower than a 
threshold velocity. Otherwise, the interplay between the geometry and causality 
leads to different scaling laws for the density of vortices as a function of the quench 
rate, as we also illustrate for the case of a toroidal trapping potential. 
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1. Introduction 

Non-equilibrium phase transitions generally lead to phases with limited long-range 
order. When a system is quenched through a critical point of a second-order phase 
transition both the correlation length £ and relaxation time r diverge. At the freeze- 
out time, t, the relaxation time of the system equals the time scale of the quench, the 
dynamics essentially freezes (impulse stage) , and there is a breakdown of adiabaticity. 
The paradigmatic Kibble-Zurek theory predicts that the average size of the domains 
in the low-symmetry phase is given by the correlation length at the freeze-out time 
[1,2] 

£:=£(t). (1) 

Ultimately, the freeze-out time depends on the quench rate 1/tq, which leads to a 
scaling law of the density of defects as a function of the rate at which the transition 
is crossed, 

v = c D ~-r%, (2) 

where D is the dimension of the domains being considered and a < 0. It is well 
known that dissipation as well as other non-universal mechanisms for defect losses 
(such as annihilation by scattering of defects with opposite topological charges) 
might lead to deviations from the Kibble-Zurek scaling whenever they become 
dominant. Nonetheless, this scenario is supported by different numerical studies [3], 
and experiments aimed at the confirmation of this prediction have been carried out in 
a variety of systems [4] , see [5] for a recent review. 

Recent experiments with pancake shaped Bose-Einstein condensates have 
reported the spontaneous nucleation of vortices during condensation [6]. Starting 
with a thermal gas in an oblate harmonic trap a linear quench in the temperature was 
applied to induce condensation. Such scenario can be naturally discussed in the light 
of the Kibble-Zurek mechanism [1, 2, 7]. A crucial feature in the experiments is the 
inhomogeneous character of the system arising from the external trapping potential. 
As a consequence, the transition does not occur simultaneously in the entire system 
and the homogeneous Kibble-Zurek mechanism (HKZM) described above breaks down, 
making it necessary to extend it to scenarios where the nucleation of defects is governed 
by causality [8, 9, 10, 11, 12, 13, 14]. To date, there is no experimental evidence 
supporting this extension that we shall refer to as the Inhomogeneous Kibble-Zurek 
mechanism (IKZM), and whose main prediction is the existence of two regimes: a) 
one in which the phase transition is crossed adiabatically with complete suppression of 
nucleation of defects b) other, characterized by vortex nucleation, where the density of 
defects after the quench obeys a scaling law with the quenching rate, different from that 
in the homogeneous mechanism described above and governed by the inhomogeneities 
in the system. Due to the high-control of the trapping potential which induces the 
inhomogeneous density profile in a trapped cloud, the nucleation of vortices during 
Bose-Einstein condensation stands out as an ideal scenario to test the predictions of 
the IKZM and it is highly desirable to extend the results of the experiments in [6] to 
such aim. Here, we analyze theoretically in this experimental setup the benchmarks 
of the IKZM, which are key for its verification. 
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2. The Homogeneous Kibble-Zurek mechanism 

We start recalling the results of the HKZM for a uniform thermal gas [15]. Consider a 
uniform quench of the temperature T(t) across the critical value of condensation T c . 
For a symmetric linear quench between the initial Tj = T c (0)+S and final Tf = T c (0)—S 
temperatures, it follows that 

T(t)=T i -t^ f -, 

T 

= T c (0)(l-±), (3) 



TQ 



where t — t — r/2 and 



The reduced temperature 

T c (0) - T(t) 

e(t) = r c (o) (5) 

governs the divergence of both the correlation length 

m - WW (6) 

and the relaxation time 

r(t) = , ■ T ° (7) 

as the system approaches the critical point (e(t) = 0). Here, {z,^} are the critical 
exponents determined by the universality class to which the system belongs to. The 
instant t in which the relaxation time equals the time remaining to the transition, 

r(t) = ^=:t, (8) 
e t 

is the freeze-out time which fixes the area £ 2 = £(t) 2 of the spots where the phase 
of the condensate is picked homogeneously. It has been theoretically shown [16] and 
experimentally demonstrated [17] that merging independent condensates with uniform 
random phases can lead to the nucleation of vortices. The geometrical configuration in 
the merging process determines the yield according to the geodesic principle [18]. The 
efficiency of this process can be captured by a constant / independent of the critical 
exponents of the system. Moreover, since the system is homogeneous, so it is the 
transition, and defects might nucleate everywhere in the system. It follows that the 
density of vortices that spontaneously nucleate under such a quench can be estimated 
as the inverse of the square of the correlation length at the freeze out time, 



Ohkzm - — - ( ^F77^ ) • (9) 



11/ t 28 

Experiments on the Bose-Einstein condensation of a 3D thermal cloud [19] have 
reported a critical exponents in agreement with the static 3D XY universality class, 
for which the best theoretical estimate to date is v — 0.6717(1) [20]. For our purposes 
the approximation v ~ 2/3 will suffice. The dynamic critical exponents is expected to 
be z = 3/2 as in the superfluid transition in 4 He, the model F in the classification of 
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Hohenberg and Halperin [21], up to possible small deviations discussed in [22]. This 
leads to a dependence 

2/3 1/2 

£>hkzm ~ t q i while £>hkzm ~ t q (10) 

follows from the mean- field values v = 1/2, z = 2. Finite-size effects might pave the 
way to an adiabatic transition whenever the correlation length at the freeze-out time 
surpasses the size of the system. Other than that, nucleation of vortices will take place 
no matter how slowly the temperature is quenched. 



3. The Inhomogeneous Kibble-Zurek mechanism 

In the following we focus on the role of the inhomogeneities arising as a result of the 
external trapping potential. We shall see that its presence brings two new ingredients, 
a local critical temperature and a local quench rate, changing the power-law for the 
density of defects as a function of the quench rate. 

Let us consider a thermal gas confined in an 3D oblate harmonic trap, isotropic 
in the radial direction, along which the density distribution exhibits a Gaussian profile 

„J 2 

of the form n(r) = n^e 2k B TT . And let us focus on the nucleation of vortices on the 
equatorial plane. Due to the inhomogeneous density, the critical temperature acquires 
a dependence on the radial position [11], 

T c {r) =T c (0)e-^r r2 . (11) 

Here, T c (0) is the critical temperature for the homogeneous system [23, 24]. For 
compactness, we shall introduce the thermal length A = y/iksT /2mLu 2 . Due to the 
spatial dependence of T c (r) , as the temperature is quenched different parts of the 
system undergo condensation at different times. It turns out to be convenient to 
introduce the spatially-dependent reduced temperature 

to identify a front in the system crossing the transition at a given position and time 
(rp,tp) satisfying the condition 

e(r F ,t F )=0. (13) 
It follows that — = 1 — 54ht, which in turn, allows us to rewrite the relative 

TQ T c (0)' ' 

temperature as e(r, t) = , with a radial-dependent quench time 

n,(D-^. (14) 

As a result, the relaxation time f = f (r) and correlation length £ — £(r) at the freeze- 
out time acquire a local dependence. The process resembles closely the formation of 
solitons in a ID BEC [11]. The nucleation of topological defects in an inhomogeneous 
transition is governed by causality [8, 9, 10, 11, 12, 13, 14]. Indeed, when the front of 
the transition moves faster than the characteristic velocity of a perturbation, defects 
nucleate, while otherwise the choice of the order parameter in the broken symmetry 
phase is done homogeneously along the system. 

Approaching the transition from the high symmetry phase, the thermal gas starts 
to condense in the center of the cloud where the density and critical temperature are 
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higher. The condensate, where the U(\) symmetry is broken, grows then radially, the 
velocity of the front being 

-l 



v F = 



Tc(0) 

T Q 
A 2 



dT c (r) 



dr 



, - (15) 
\r\r Q (r) 

where we have disregarded corrections due to the time derivative of the critical 
temperature T c (r) — T c {r, t) arising from the quench T(t) which affects the density 
profile. This is reasonable whenever the amplitude of the quench is small 5 <C 1^(0). 
The existence of this front and a local temperature play a crucial role in the following, 
and requires a local thermalization fast compared to the time scale of the quench. For 
defects to nucleate, the front of the transition is to move faster than any perturbation. 
The characteristic velocity of a perturbation when the dynamics stops being adiabatic, 
this is, at the freeze-out time, can be upper-bounded by the ratio of the freeze-out 
correlation length over the relaxation time, 

i & ( t 



6 = t = — ( — 7-T J • (16) 
T r o \TQ(r)J 

While within the HKZM the appearance of vortices is always expected (up to finite-size 
effects of the system), the trapping potential paves the way for a perfectly adiabatic 
condensation, where the spontaneous nucleation is suppressed as long as 

v F <v. (17) 

Conversely, the condition for appearance of defects, v F > v, leads to a threshold value 
of the radius of the cloud delimiting the area within which defects might nucleate, 



< 



Co \TQ(r) 



A 2 ( t \ l +<"> 



<=-(^) exp(^^l). (18) 



AVV\^ ( r 2 1 
Co {r Q J CXP UAM + ^ 

Within a reduced area of the condensate <S*, formation of vortices occurs as in the 
HKZM, setting the average size £ 2 of the spots with uniform local phase. Up to a 
numerical factor /, the total number of topological defects cab be estimated as follows 

r 2nr 

Mkzm - / dr TI7^ (I 9 ) 

J{r\v F >v} f£(r) 2 

and simply by A/ikzm = jp if one neglects the radial dependence of £. 

The transcendental inequality (18) can be approximated around the center of the 
cloud (f <C A) by setting tq(t-) = tq(0) = tq, and consistently £(r) = £(0). Solving 
the associate equality one finds the effective radius f = r* , in terms of which 5* = nr 2 
so that 

(0) _^f^_\ 2J ^ 
A/ikzm- ^ [ tTc{q) ) ■ (20) 

We note that for the mean-field critical exponents, z = 2 and v = 1/2, the exponent 
governing the scaling in Eq. (20) is given by 2 ^+^ — 2, which is four times that 
predicted by the HKZM. For z — 3/2 and v = 2/3, the IKZM power-law exponent 




Figure 1. Ratio of the velocity of the front vp crossing the critical point for 
condensation, and the characteristic speed of a perturbation v at the freeze- 
out time, along the radial axis of a pancake-shaped atomic thermal cloud, 
for the values A = 0.1, 0.4, l/\/e, 0.8 (from bottom to top). The coefficient 

= \J i+tl x+uz allows one to identify the homogeneous (A > 1/y/e) and 
inhomogeneous (A < 1/y/e) scenarios. 



becomes 7/3, which is as well nearly four times that of the homogeneous counterpart, 
and constitutes an experimentally accessible benchmark of the IKZM. 

Nonetheless, note that the absolute number of defects is reduced by a factor 



/ A 2 \ * / 
( A- \ (T 



2(1 + ./) 

Mkzm/A/hkzm ~ {n/r M f ~ ( 7^ ) ( — ) ^ , (21) 

\t,or M J \tqJ 

where tm is the radius of the cloud. Hence, it is the dependence on the quench rate 
and the critical exponents inherited by the effective size of the cloud iS* and £ that is 
responsible for the new power-law governing the density of defects. 

For the general solution of Eq. (18), it is convenient to introduce the variable 



to find £ < Ae 2 , with 



1+" 



Co V 1 + vz \tq 

Figure 1 shows the ratio vf/v along the radial coordinate of the cloud for different 
values of this parameter. In the central region where the density reaches its maximum, 
and hence it is approximately uniform, this ratio is always larger than unity. 

One can distinguish two different regimes as a function of the value of A with 
respect to the critical value A c — l/\/e: (i) For A > A c , vf is everywhere along the 
sample larger than v, so the homogeneous KZM applies. The density of vortices is then 
given by Eq. (9). (ii) For A < A c there exist two solutions {£* = C{ r *) < C* = C(K)} 
defining two disjoint concentric disks with support on the interval T :— [0, r*] U [r^, ry[ ] 
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Figure 2. 

Reduction due to causality of the effective size of the cloud for nucleation 
of vortices, a) Critical values r» and r' t of the radial coordinate of the cloud 
determining the regions where defects might nucleate as a function of A < Ac- b) 
Effective area as a function of the parameter A within the IKZM for a pancake 
atomic cloud as well as for a cloud in a toroidal trap, which will be discussed in 
Section 4. 



along the radial direction, where rjvf is the effective radius of the cloud. As a result, 
see Fig. 2 , the area 5* (A) = ir [rl + {r 2 M — r' 2 )] where vortices might nucleate is reduced 
with respect the total area of the cloud S = irr 2 M , and so it is the corresponding density 
of defects (note that case (i) corresponds to r* = r^, so the effective area equals 

the total area S). Nonetheless, for A < Am = C( r M)e~ " , r* > tm so the outer 
disk can be ignored, and nucleation can be considered to take place only in the center 
of the cloud (r < r*). This is the regime where Eq. (20) holds, up to finite size effects 
which might lead to an adiabatic transition and a breakdown of the scaling whenever 
£ > r* . These different cases follow from the spatial distribution of the regions where 
the ratio Vf/£>lasa function of A, exhibiting a transition from case (i) to (ii). For 
instance, if in a given experiment only the quenching time is varied, one can identified 
a critical value 

(24) 

around which the scaling changes between those predicted by HKZM and IKZM. 

To appreciate how the scaling is modified, we introduce the dimensionless cloud 
width and quench rate 

e = A/£ , T = t /tq, (25) 

and rewrite the parameter A as 

A= x l\±^m^. (26) 



tq(A c ) = T 



Co 



1 + vz 



y/eA V 1 + 2i 
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Figure 3. Density of vortices as a function of the dimensionless width of the cloud 
and quench rate T in an inhomogeneous phase transition under a linear quench 
of the temperature a) in a radially symmetric harmonic trap b) in a toroidal 
trap. Here v = 2/3, z = 3/2, but the same qualitative behaviour is observed for 
v = 1/2, z = 2. For small values of (0, T) the density of defects exhibits the 
IKZM power-law scaling of Eq. (20) in a harmonic confinement, and that of Eq. 
(30) in a toroidal trap. As (0, T) are increased, within the IKZM the power-law 
dependence breaks down and exhibits a more complicated dependence. For large 
enough values of 0, and T (such that A > -4 C ), nucleation is possible in the whole 
cloud, the confinement does not play a role any more, an the HKZM scaling in 
Eq. (9) describes the dependence of the density of vortices in both type of traps. 

Equating both sides of the inequality Eq. (18), one can find the values of r*, r£ that 
determine the effective area <S*(9, T) = S*(A)/A 2 where vf/v > 1. This leads to the 
following equations, 

C.(6, T) = J±±^or^e< 2 *( & - r » 2 , (27) 

M K zM(e,r) = ^i±^e 2 T^ f dcce^, (28) 

/ 1 + 2^ J{c\v F >v} 

which take into account the causality argument and the local dependence of £. Though 
this expression lacks a simple power-law scaling with the temperature quench and 
sweeping rate through the transition, it is still universal, in the sense that such 
dependencies are still governed by the critical exponents associated with universality 
class to which the transition belongs. Figure 3a shows the density of topological 
defects T> = Af/S as a function the dimensionless temperature and rate. As a result 
of the inhomogeneity the dependence of the density of vortices on the cooling rate 
becomes much stronger than in the homogeneous case. Moreover, the threshold Am 
that arises from the finite-size of the cloud, is responsible for an abrupt jump on 
density of vortices. For the sake of illustration in Fig. 3, we consider the case in which 
n(rM)/n(0) = O.f, but the results are qualitatively the same for different values of 
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In the experiments reported in [6], a fitting to the measured temperatures to 
(1 — t/rg) leads to an estimate tq ~ 5 s, while Tq ~ 0.1 s is given by the scattering 
time of atoms. Moreover, the de Broglie wavelength is £ ~ 1-6 Atrn and A ~ 65 /um, 
so that A <~ 1.7, suggesting a homogeneous scenario in the proximity of the boundary 
A c - Using the explicit form of A as a guide, we note that an experimental verification 
of the IKZM, and in particular of the scaling in Eq. (20), will be favoured by a slow 
quench, T = tq/tq < 1, and a tight trapping potential so as to reduce = A/£o- 

We close this section noticing that our discussion is limited to quenches where 
the correlation length at the freeze out time is much smaller than the trap size. 
Otherwise, the finite-size scaling modifies the power-law governing the divergence of 
the correlation length in the neighbourhood of the critical point, as discussed in [25]. 



4. Mixing the mechanisms: Bose-Einstein condensation in a toroidal trap 

The Tucson group also studied the spontaneous nucleation of vortices in a toroidal 
trap, that we describe in the following. Such geometry induces a new type of transition 
with a mixed homogeneous-inhomogeneous character. To illustrate this, consider 
a squashed toroidal trap of radius r c , and transverse trapping frequency u r . The 
corresponding equilibrium density profile of a non-interacting thermal gas has the form 
n(9, r) = n exp — (r — r c ) 2 /2A 2 , and as a result the critical temperature becomes 



T c (9,r) =T c exp 



(r-r c ) 2 ) 
2A 2 



(29) 



The upshot is that the transition remains homogeneous as a function of 9 but it is 
governed by causality in the radial direction where it is still inhomogeneous. The 
relative coordinate h = r — r c plays the role of r in the preceding section. As in 
the harmonic trap, a complete suppression of vortex nucleation can arise not only 
from finite-size effects, but also whenever vf < v. As for the pancake condensate, A c 
governs the transition between the homogeneous and inhomogeneous scenario in the 
radial direction. For A < A c the effective size of the cloud is reduced as shown in Fig. 
2b, and the IKZM scaling of defects is modified as follows 

(0) _ 47TTcfe. _ 4^r c A 2 / t 2S \ ^ 
JVikzm- ^ - ^3 ^ rTc(0) J > (30) 

under the approximation mentioned above, /i* = |r* — r c \ <C A, £ = £(r c ), which 

holds for A < C(l r c — f~M|)e 5 . For the mean-field critical exponents 
~ t q ' , at variance with the HKZM power law with exponent 1/2 in Eq. 

(9) and the inhomogeneous case, see Eq. (20) where the exponent is ^yr^-^ = 2. 
Similarly, for v = 2/3, z = 3/2, the power-law exponent is 3/2 intermediate between 
that of HKZM (2/3) and IKZM (7/3). The general estimate in the inhomogeneous 
scenario (A > A c ) can be obtained by solving numerically Eq. (27) as a function of the 
dimcnsionless transverse width of the cloud & and quench rate T, and noticing that the 
effective area is in this case given by 5* (9, T) = 47r£(r c )[C*(9, T) — T) + C( r M)]- 
The result is shown in Fig. (3)b, which illustrates the change in the scaling from 
that of Eq. (9) (HKZM) to that in Eq. (30) (IKZM) as the quenching time and the 
inhomogeneity of the cloud are increased. The boundary between the different scalings 
is now smoother than for the pancake geometry, due to the mixed character of the 
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transition. We further note that for tight toroidal traps even if vortex nucleation 
is suppressed, solitons might be formed leading to the spontaneous generation of 
persistent currents as recently discussed in [26]. Nonetheless, for wider traps such 
solitons become unstable against vortex formation through the snake instability and 
related mechanisms [27]. 

5. Conclusions 

In conclusion, we have shown of the Kibble-Zurek mechanism, extended to describe 
spatially inhomogeneous systems, severely modifies with respect to the homogeneous 
case the scaling of the density of defects as a function of the quenching rate, 
as illustrated in Table 1. When the presence of a trapping potential leads to a 



Critical exponents \ Trap 


Homogeneous 


Harmonic 


Toroidal 


Arbitrary (v, z) 
Mean- field theory (v = |, z = 2) 
Experiments/F model (v = |, z = §) 


2v 
l+vz 

l 

2 

2 
3 


2{l+2v) 
1 + vz 

2 

7 

3 


1+3;/ 
l + vz 

5 
4 

3 
2 



Table 1. Power laws predicted by the Kibble-Zurek mechanism for the density 
of vortices T> as a function of the rate tq of a thermal quench through the critical 
temperature for Bosc-Einstcin condensation. The exponent a of the power-law 
T> ~ t q C " i s shown for different critical exponents (u, z) and trapping potentials. 
The scalings in the harmonic toroidal and toroidal trap are restricted to situations 
where nucleation is limited to a small fraction of the cloud (r*,h* <C 1) as 
discussed in the text. 

inhomogeneous scenario, a neat power-law scaling governs the nucleation of vortices 
only when causality limits it to a small fraction of the cloud. Otherwise, the local 
dependence of the effective quench rate is to be taken into account, leading to a 
more complicated behaviour, different from a power-law. We have further introduced 
a simple parameter A which allows one to estimate which is the relevant scenario 
for the nucleation of defects in a trapped cloud. Indeed, the standard Kibble-Zurek 
scaling is recovered for weak trapping potentials and fast quenches. We close noticing 
that while the homogeneous Kibble-Zurek mechanism has been studied in a wide 
variety of experiments [4] , the mechanism for the nucleation of topological defects in 
inhomogeneous systems lacks to date experimental evidence. Hence, we hope that 
measurements of the number of spontaneously generated vortices as a function of the 
cooling rate of an atomic cloud undergoing Bose-Einstein condensation, might soon 
change this state of affairs. 
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